#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <unistd.h>
#include <math.h>
#include <float.h>
#include"gps.h"

/*********************************************************************
	YMD_To_MJD - Compute the Modified Julian Date (MJD) from year, month, and day

	MJD = Julian Date - 2400000.5
**********************************************************************/
double YMD_To_MJD(
	int iYear, //i -year
	int iMonth, //i - month
	double dfDay //i day of month
	)
{
	int J,M,I,K;
	double MJD;

	J=iYear;
	M=iMonth;

	if(M>2){
		I=J/100;
		K=2-I+(I/4);
		MJD=(365.25*J-fmod(365.25*J,1.0))-679006.0;
		MJD+=(int)(30.6001*(M+1))+dfDay+K;
	}else{
		J=J-1;
		M=M+12;
		I=J/100;
		K=2-I+(I/4);
		MJD=(365.25*J-fmod(365.25*J,1.0))-679006.0;
		MJD+=(int)(30.6001*(M+1))+dfDay+K;
	}
	return MJD;
}

/*********************************************************************
	MJD_To_YMD - Compute year, month, day of month from the Modified Julian Date (MJD)  

	MJD = Julian Date - 2400000.5
**********************************************************************/
void MJD_To_YMD(
	double dfMJD,	//i - Modified Julian Date
	int* piYear,	//o - year
	int* piMonth,	//o - month
	double* pfDay	//o - day of month
	)
{
	double T1,T2,T3,T4;
	int IH,IH1,IH2;

	T1=1.0+dfMJD -fmod(dfMJD,1.0)+2400000.0;
	T4=fmod(dfMJD,1.0);
	IH=(int)((T1-1867216.25)/36524.25);
	T2=T1+1+IH-(IH/4.0);
	T3=T2-1720995.0;
	IH1=(int)((T3-122.1)/365.25);
	T1=365.25*IH1-fmod(365.25*IH1,1.0);
	IH2=(int)((T3-T1)/30.6001);
	*pfDay=(T3-T1-(int)(30.6001*IH2))+T4;
	*piMonth=(int)(IH2-1.0);

	if(IH2>13){
		*piMonth=IH2-13;
	}
	
	*piYear=IH1;

	if(*piMonth<=2){
		*piYear=*piYear+1;
	}
}

void test_mjd1(void)
{
	int iYear=2019;
	int iMonth=4;
	double dfDay=29;
	double mjd;

	int hour;
	int minute;
	int second;

	hour=10;
	minute=45;
	second=48;
	printf("iYear=%d,iMonth=%d,dfDay=%.40g,hour=%d,minute=%d,second=%d\n",
			iYear,iMonth,dfDay,hour,minute,second);

	dfDay+=(double)second/3600.0/24.0;
	dfDay+=(double)minute/60.0/24.0;
	dfDay+=(double)hour/24.0;

	mjd=YMD_To_MJD( iYear, iMonth, dfDay);
	printf("mjd=%.40g\n",mjd);

	iYear=0;
	iMonth=0;
	dfDay=0;
	hour=0;
	minute=0;
	second=0;
			
	MJD_To_YMD(mjd,&iYear,&iMonth,&dfDay);
	mjd=dfDay-floor(dfDay);
	dfDay=floor(dfDay);
	mjd*=24;
	hour=floor(mjd);
	mjd-=hour;
	mjd*=60;
	minute=floor(mjd);
	mjd-=minute;
	mjd*=60;
	second=round(mjd);
	
	printf("iYear=%d,iMonth=%d,dfDay=%.40g,hour=%d,minute=%d,second=%d\n",
			iYear,iMonth,dfDay,hour,minute,second);	
}


void test_bd_minus_gps(void)
{
	int iYear;
	int iMonth;
	double dfDay;

	int hour;
	int minute;
	int second;

	double mjd_bd;
	double mjd_gps;
	double mjd;

	iYear=2006;
	iMonth=1;
	dfDay=1;
	hour=0;
	minute=0;
	second=0;


	dfDay+=(double)second/3600.0/24.0;
	dfDay+=(double)minute/60.0/24.0;
	dfDay+=(double)hour/24.0;

	mjd_bd=YMD_To_MJD( iYear, iMonth, dfDay);

	iYear=1980;
	iMonth=1;
	dfDay=6;
	hour=0;
	minute=0;
	second=0;


	dfDay+=(double)second/3600.0/24.0;
	dfDay+=(double)minute/60.0/24.0;
	dfDay+=(double)hour/24.0;

	mjd_gps=YMD_To_MJD( iYear, iMonth, dfDay);

	mjd=mjd_bd-mjd_gps;
	printf("mjd=%.40g\n",mjd);

	iYear=0;
	iMonth=0;
	dfDay=0;
	hour=0;
	minute=0;
	second=0;
			
//	MJD_To_YMD(mjd,&iYear,&iMonth,&dfDay);
	mjd=dfDay-floor(dfDay);
	dfDay=floor(dfDay);
	mjd*=24;
	hour=floor(mjd);
	mjd-=hour;
	mjd*=60;
	minute=floor(mjd);
	mjd-=minute;
	mjd*=60;
	second=round(mjd);
	
	printf("iYear=%d,iMonth=%d,dfDay=%.40g,hour=%d,minute=%d,second=%d\n",
			iYear,iMonth,dfDay,hour,minute,second);	
}

int main(int argc, char *argv[])
{
	test_bd_minus_gps();

	return 0;
}

